// hnum_cblas.h
#pragma once
#ifndef __HNUM_CBLAS_H__
#define __HNUM_CBLAS_H__

#include "hnumdef.h"

#include <math.h>

#include "hnum_slu_mt_util.h"
#include "hnum_slu_scomplex.h"
#include "hnum_slu_dcomplex.h"
#include "hnum_supermatrix.h"
#include "hnum_pxgstrf_synch.h"

namespace harlinn
{
    namespace numerics
    {
        namespace SuperLU
        {

            typedef int integer;
            typedef int logical;

            typedef char *address;
            typedef float real;
            typedef double doublereal;
            typedef short int shortlogical;
            typedef char logical1;
            typedef char integer1;

            inline int lsame_(const char *first, const char *second)
            {   
                if(*first == *second)
                {
                    return true;
                }
                return toupper(*first) == toupper(*second);
            }

            inline int lsame_(const char first, const char second)
            {   
                if(first == second)
                {
                    return true;
                }
                return toupper(first) == toupper(second);
            }

            int xerbla_(char *, int *);


            // BLAS

            // constant times a vector plus a vector
            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int caxpy_(integer *n, complex *ca, complex *cx, integer * incx, complex *cy, integer *incy);

            // constant times a vector plus a vector
            // In use if USE_VENDOR_BLAS is *not* defined
            inline int caxpy_(integer n, const complex *ca, const complex *cx, complex *cy)
            {
                if ((n <= 0)||( fabs(ca->r) + fabs(ca->i) == 0.f))
                {
	                return 0;
                }
                for (int i = 0; i < n; i++) 
                {
	                cy[i].r = cy[i].r + (ca->r * cx[i].r - ca->i * cx[i].i);
                    cy[i].i = cy[i].i + (ca->r * cx[i].i + ca->i * cx[i].r);
                }
                return 0;
            }

            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int ccopy_(integer *n, complex *cx, integer *incx, complex *cy, integer *incy);

            // In use if USE_VENDOR_BLAS is *not* defined
            inline int ccopy_(integer n, complex *cx, complex *cy)
            {
                memcpy(cy,cx,n * sizeof(complex));
                return 0;
            } 


            // forms the dot product of two vectors, conjugating the first vector.
            // Not used
            HNUM_EXPORT void cdotc_(complex * ret_val, integer *n, complex *cx, integer *incx, complex *cy, integer *incy);

            // forms the dot product of two vectors, conjugating the first vector.
            // Not used
            inline void cdotc_(complex * ret_val, integer n, const complex *cx, const complex *cy)
            {
                complex ctemp = {0,0};
                if (n <= 0) 
                {
                    *ret_val = ctemp;
	                return;
                }
                for (int i = 0; i < n; i++) 
                {
	                ctemp.r += cx[i].r * cy[i].r + cx[i].i * cy[i].i;
                    ctemp.i += cx[i].r * cy[i].i - cx[i].i * cy[i].r;
                }
                *ret_val = ctemp;
            }


            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int cgemv_(char *trans, integer *m, integer *n, complex * alpha, complex *a, integer *lda, complex *x, integer *incx, complex *beta, complex *y, integer *incy);
            // Not used
            HNUM_EXPORT int cgerc_(integer *m, integer *n, complex *alpha, complex *x, integer *incx, complex *y, integer *incy, complex *a, integer *lda);

            // Not used
            HNUM_EXPORT int chemv_(char *uplo, integer *n, complex *alpha, complex *a, integer *lda, complex *x, integer *incx, complex *beta, complex *y,integer *incy);

            // Not used
            HNUM_EXPORT int cher2_(char *uplo, integer *n, complex *alpha, complex *x, integer *incx, complex *y, integer *incy, complex *a, integer *lda);

            // Not used
            HNUM_EXPORT int cscal_(integer *n, complex *ca, complex *cx, integer *incx);

            // In use
            HNUM_EXPORT int ctrsv_(char *uplo, char *trans, char *diag, integer *n, complex *a, integer *lda, complex *x, integer *incx);


            // In use
            HNUM_EXPORT doublereal dasum_(integer *n, doublereal *dx, integer *incx);

            // In use
            HNUM_EXPORT int daxpy_(integer *n, doublereal *da, doublereal *dx, integer *incx, doublereal *dy, integer *incy);

            // In use
            HNUM_EXPORT doublereal dcabs1_(doublecomplex *z);

            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int dcopy_(integer *n, doublereal *dx, integer *incx, doublereal *dy, integer *incy);

            // In use if USE_VENDOR_BLAS is *not* defined
            inline int dcopy_(integer n, const doublereal *dx, doublereal *dy )
            {
                memcpy(dy,dx,(n)*sizeof(doublereal));
                return 0;
            }

            // Not in use
            HNUM_EXPORT doublereal ddot_(integer *n, doublereal *dx, integer *incx, doublereal *dy, integer *incy);

            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int dgemv_(char *trans, integer *m, integer *n, doublereal *alpha, doublereal *a, integer *lda, doublereal *x, integer *incx, doublereal *beta, doublereal *y, integer *incy);

            // Not in use
            HNUM_EXPORT int dger_(integer *m, integer *n, doublereal *alpha, doublereal *x, integer *incx, doublereal *y, integer *incy, doublereal *a, integer *lda);

            // Not in use
            HNUM_EXPORT doublereal dnrm2_(integer *n, doublereal *x, integer *incx);

            // Not in use
            HNUM_EXPORT int drot_(integer *n, doublereal *dx, integer *incx, doublereal *dy, integer *incy, doublereal *c, doublereal *s);

            // Not in use
            HNUM_EXPORT int dscal_(integer *n, doublereal *da, doublereal *dx, integer *incx);

            // Not in use
            HNUM_EXPORT int dsymv_(char *uplo, integer *n, doublereal *alpha, doublereal *a, integer *lda, doublereal *x, integer *incx, doublereal *beta, doublereal *y, integer *incy);

            // Not in use
            HNUM_EXPORT int dsyr2_(char *uplo, integer *n, doublereal *alpha, doublereal *x, integer *incx, doublereal *y, integer *incy, doublereal *a, integer *lda);

            // In use
            HNUM_EXPORT int dtrsv_(char *uplo, char *trans, char *diag, integer *n, doublereal *a, integer *lda, doublereal *x, integer *incx);

            // Not in use
            HNUM_EXPORT doublereal dzasum_(integer *n, doublecomplex *zx, integer *incx);

            // Not in use
            HNUM_EXPORT doublereal dznrm2_(integer *n, doublecomplex *x, integer *incx);

            // Not in use
            HNUM_EXPORT integer icamax_(integer *n, complex *cx, integer *incx);

            // In use
            HNUM_EXPORT integer idamax_(integer *n, doublereal *dx, integer *incx);

            // In use
            HNUM_EXPORT integer isamax_(integer *n, real *sx, integer *incx);

            // Not in use
            HNUM_EXPORT integer izamax_(integer *n, doublecomplex *zx, integer *incx);

            // In use
            HNUM_EXPORT real sasum_(integer *n, real *sx, integer *incx);

            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int saxpy_(integer *n, real *sa, real *sx, integer *incx, real *sy, integer *incy);

            // constant times a vector plus a vector.   
            // In use if USE_VENDOR_BLAS is *not* defined
            inline int saxpy_(integer n, const real *sa, const real *sx, real *sy)
            {
                if (n <= 0) 
                {
	                return 0;
                }
                if (*sa == 0.f) 
                {
	                return 0;
                }
                int m = n % 4;
                int i = 0;
                if (m != 0) 
                {
                    for (;i < m; i++) 
                    {
	                    sy[i] += *sa * sx[i];
                    }
                    if (n < 4) 
                    {
	                    return 0;
                    }
                }
                for (; i < n; i += 4) 
                {
	                sy[i]     += *sa * sx[i];
	                sy[i + 1] += *sa * sx[i + 1];
	                sy[i + 2] += *sa * sx[i + 2];
	                sy[i + 3] += *sa * sx[i + 3];
                }
                return 0;
            }

            // Not in use
            HNUM_EXPORT real scasum_(integer *n, complex *cx, integer *incx);

            // Not in use
            HNUM_EXPORT real scnrm2_(integer *n, complex *x, integer *incx);

            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int scopy_(integer *n, real *sx, integer *incx, real *sy, integer *incy);

            // In use if USE_VENDOR_BLAS is *not* defined
            inline int scopy_(integer n, const real *sx, real *sy )
            {
                memcpy(sy,sx,n*sizeof(real));
                return 0;
            }

            // Not in use
            HNUM_EXPORT real sdot_(integer *n, real *sx, integer *incx, real *sy, integer *incy);

            // In use
            HNUM_EXPORT int sgemv_(char *trans, integer *m, integer *n, real *alpha, real *a, integer *lda, real *x, integer *incx, real *beta, real *y, integer *incy);

            // Not in use
            HNUM_EXPORT int sger_(integer *m, integer *n, real *alpha, real *x, integer *incx, real *y, integer *incy, real *a, integer *lda);

            // Not in use
            HNUM_EXPORT real snrm2_(integer *n, real *x, integer *incx);

            // Not in use
            HNUM_EXPORT int srot_(integer *n, real *sx, integer *incx, real *sy, integer *incy, real *c, real *s);

            // Not in use
            HNUM_EXPORT int sscal_(integer *n, real *sa, real *sx, integer *incx);

            // Not in use
            HNUM_EXPORT int ssymv_(char *uplo, integer *n, real *alpha, real *a, integer *lda, real *x, integer *incx, real *beta, real *y, integer *incy);
            // Not in use
            HNUM_EXPORT int ssyr2_(char *uplo, integer *n, real *alpha, real *x, integer *incx, real *y, integer *incy, real *a, integer *lda);

            // In use
            HNUM_EXPORT int strsv_(char *uplo, char *trans, char *diag, integer *n, real *a, integer *lda, real *x, integer *incx);

            // Not in use
            HNUM_EXPORT int zaxpy_(integer *n, doublecomplex *za, doublecomplex *zx, integer *incx, doublecomplex *zy, integer *incy);

            // Not in use
            inline int zaxpy_(integer n, doublecomplex *za, doublecomplex *zx, doublecomplex *zy)
            {
                if ((n <= 0)||( fabs(za->r) + fabs(za->i) == 0.f))
                {
	                return 0;
                }
                for (int i = 0; i < n; i++) 
                {
	                zy[i].r = zy[i].r + (za->r * zx[i].r - za->i * zx[i].i);
                    zy[i].i = zy[i].i + (za->r * zx[i].i + za->i * zx[i].r);
                }
                return 0;
            }

            // In use if USE_VENDOR_BLAS is defined
            HNUM_EXPORT int zcopy_(integer *n, doublecomplex *zx, integer *incx, doublecomplex *zy, integer *incy);

            // In use if USE_VENDOR_BLAS is *not* defined
            inline int zcopy_(integer n, const doublecomplex *zx, doublecomplex *zy )
            {
                memcpy(zy,zx,n*sizeof(doublecomplex));
                return 0;
            }

            // Not in use
            HNUM_EXPORT void zdotc_(doublecomplex * ret_val, integer *n, doublecomplex *zx, integer *incx, doublecomplex *zy, integer *incy);

            // Not in use
            inline void zdotc_(doublecomplex * ret_val, integer n, const doublecomplex *zx, const doublecomplex *zy)
            {
                doublecomplex ztemp = {0,0};
                if (n <= 0) 
                {
                    *ret_val = ztemp;
	                return;
                }
                for (int i = 0; i < n; i++) 
                {
	                ztemp.r += zx[i].r * zy[i].r + zx[i].i * zy[i].i;
                    ztemp.i += zx[i].r * zy[i].i - zx[i].i * zy[i].r;
                }
                *ret_val = ztemp;
            }

            // In use
            HNUM_EXPORT int zgemv_(char *trans, integer *m, integer *n, doublecomplex *alpha, doublecomplex *a, integer *lda, doublecomplex *x, integer *incx, doublecomplex *beta, doublecomplex *y, integer *incy);

            // Not in use
            HNUM_EXPORT int zgerc_(integer *m, integer *n, doublecomplex *alpha, doublecomplex *x, integer *incx, doublecomplex *y, integer *incy, doublecomplex *a, integer *lda);
            // Not in use
            HNUM_EXPORT int zhemv_(char *uplo, integer *n, doublecomplex *alpha, doublecomplex *a, integer *lda, doublecomplex *x, integer *incx, doublecomplex *beta, doublecomplex *y, integer *incy);

            // Not in use
            HNUM_EXPORT int zher2_(char *uplo, integer *n, doublecomplex *alpha, doublecomplex *x, integer *incx, doublecomplex *y, integer *incy, doublecomplex *a, integer *lda);

            // Not in use
            HNUM_EXPORT int zscal_(integer *n, doublecomplex *za, doublecomplex *zx, integer *incx);

            // In use
            HNUM_EXPORT int ztrsv_(char *uplo, char *trans, char *diag, integer *n, doublecomplex *a, integer *lda, doublecomplex *x, integer *incx);

            HNUM_EXPORT void clsolve ( int ldm, int ncol, complex *M, complex *rhs );
            HNUM_EXPORT void cusolve (int ldm,int ncol, complex *M,complex *rhs);
            HNUM_EXPORT void cmatvec (int, int, int, complex *, complex *, complex *);

            HNUM_EXPORT void dlsolve ( int ldm, int ncol, double *M, double *rhs );
            HNUM_EXPORT void dusolve (int ldm,int ncol, double *M,double *rhs );
            HNUM_EXPORT void dmatvec (int ldm,int nrow,int ncol,double *M,double *vec,double *Mxvec);
            HNUM_EXPORT void dmatvec2 (int lda,int m,int n,double *A,double *x0,double *x1,double *y0,double *y1);

            HNUM_EXPORT void slsolve ( int ldm, int ncol, float *M, float *rhs );
            HNUM_EXPORT void susolve ( int ldm, int ncol, float *M, float *rhs );
            HNUM_EXPORT void smatvec (int ldm,int nrow,int ncol,float *M,float *vec,float *Mxvec);
            HNUM_EXPORT void smatvec2 (int lda,int m,int n,float *A,float *x0,float *x1,float *y0,float *y1);

            HNUM_EXPORT void zlsolve ( int ldm, int ncol, doublecomplex *M, doublecomplex *rhs );
            HNUM_EXPORT void zusolve ( int ldm,	int ncol, doublecomplex *M,	doublecomplex *rhs );
            HNUM_EXPORT void zmatvec ( int ldm,	int nrow,  int ncol, doublecomplex *M,	doublecomplex *vec,	doublecomplex *Mxvec );


            


            // LAPACK 

            HNUM_EXPORT double dzsum1_(int *n, doublecomplex *cx, int *incx);
            HNUM_EXPORT double scsum1_(int *n, complex *cx, int *incx);

            HNUM_EXPORT double dlaran_(int *iseed);
            //
            HNUM_EXPORT void heap_relax_snode (const int n,superlumt_options_t *superlumt_options, pxgstrf_relax_t *pxgstrf_relax );


            // 
            HNUM_EXPORT int genmmd_(int *neqns, int *xadj, int *adjncy, 
	                        int *invp, int *perm, int *delta, int *dhead, 
	                        int *qsize, int *llist, int *marker, int *maxint, 
	                        int *nofsub);

            HNUM_EXPORT int clacon_(int *n, complex *v, complex *x, float *est, int *kase);
            HNUM_EXPORT int dlacon_(int *n, double *v, double *x, int *isgn, double *est, int *kase);
            HNUM_EXPORT int slacon_(int *n, float *v, float *x, int *isgn, float *est, int *kase);
            HNUM_EXPORT int zlacon_(int *n, doublecomplex *v, doublecomplex *x, double *est, int *kase);

            HNUM_EXPORT int icmax1_(int *n, complex *cx, int *incx);
            HNUM_EXPORT int izmax1_(int *n, doublecomplex *cx, int *incx);
            

            
        };
    };
};

#endif //__HNUM_CBLAS_H__